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1 Introduction 



^ When particles are suspended in an homogeneous isotropic fluid, 

O the viscosity of the resulting complex fluid is increased. In the 

I — I case of dilute suspensions, the increase in viscosity as a function 

^—1 of the volume fraction was firstly determined for spherically 

t5«. shaped particles by Einstein in 1911, [Ij. Some years later, ex- 

tensions of Einstein's work appeared and new formulas for the 
viscosity as a function of the volume fraction for sohd elhpsoidal 
O particles and emulsions were also derived by Jeffery (1922) and 



Taylor (1932), respectively 



•^ For sufficiently low particle concentrations, the viscosity r/ of a 

vh suspension can in general, be written as. 



r/(0) = r/o(l+M0 + A:ff02^...), (1) 

where r/o is the solvent viscosity, [77] is the low filling fraction in- 
trinsic viscosity, kn is the so called Huggins coefficient, and (j) is 
the volume fraction of the particles. The value of [77] depends on 
the particle shape. Although its calculation is difficult, there are 
few cases where [77] has been obtained analytically as in the case 



of spheres [I], ellipsoids |1], |5], long cylinders [B] and dumbbells 
consisting of two identical spheres [S] , [7] . Fortunately, accurate 
numerical approaches have been developed to calculate the in- 
trinsic viscosity of arbitrarily shaped particles [S] , p] . 

The large amount of work devoted to determine the relation be- 
tween viscosity and concentration under different physical condi- 
tions that reflect experimental protocols p!0lfTTlfT2lfT3j . is clearly 



connected with the important role that suspensions and emul- 
sions play in almost all fields of industry, medicine and biology- 
related soft-matter systems. 



Many remarkable theoretical and, more recently, numerical works 
have made central contributions to the understanding of the rhe- 
ology of concentrated suspensions (see for instance, Refs. (TH- 



I). Essentially, most of these works introduce particle correla- 
tions by taking into account hydrodynamic interactions and pro- 
vide a conceptual framework that explains how the microstruc- 
ture of the systems modifies the behavior of the viscosity as a 
function of the volume fraction. However, the quantitative de- 
scription of the problem still constitutes an open challenge for 
theoretical descriptions even in the simplest case of spherical par- 
ticles. 

Several models have been proposed in order to extend the range 
of applicability of Einstein's expression to larger volume frac- 
tions |ID]. Semi-empirical P51I5D] and effective-medium |5^ - 



models intended to extend the quantitative description of the 
experimental dependence of the effective viscosity of hard parti- 
cle suspensions at arbitrary volume fractions have been proposed 
long ago. Among them we can mention the one derived by Saito 
[H] that accounts for hydrodynamic interactions between uncor- 
related particles 



ii{(P) = m i + M T^ . (2) 



Other approaches include the use of differential effective- medium 
theories (DEMT). Krieger and Dougherty [BS] derived the follow- 
ing expression that incorporates excluded volume using a crowd- 
ing effect introduced by Mooney 



n{<p)=nJi--^] , (3) 

with 0niax the filling fraction at maximum packing. This formula 
agrees reasonably well with the experimental data particularly 
in the low concentration regime. Moreover, this relation reduces 
to the correct Einstein's equation in the limit of infinite dilu- 
tion. Very recently, Bullard and coworkers [12] also used DEMT 
techniques to derive a similar relationship for suspensions of par- 
ticles that may themselves incorporate some of the solvent either 
by solvation or by occlusion in interstitial pores. In the case of 
spheres their expression adopts the form 

r/(0)=r/o(l-i^0)-['']/^, (4) 

where the factor K considers flocculation of the particles thus 
representing the ratio of the volume of the clusters to the volume 
of the particles forming the clusters 



The original formulations of DEMT were based on Einstein's ex- 
pression, providing moderate agreement with experimental re- 
sults. It is important to stress that one common characteristic of 
these models is that they describe acceptably well the viscosity- 
concentration relation at low concentrations but fail to describe 
correctly the whole concentration regime. As we will explain later, 
this is because these theories do not incorporate appropriately the 
correlations introduced by the excluded volume effects. 

A different approach to calculate the viscosity of suspensions is 
based on the expected theoretical divergence of viscosity near the 
percolation threshold 
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T] (0) ~ r/o f 1 - ^1 , (5) 



which, according to Douglas and coworkers [32], should be inde- 
pendent of the shape of the particles in suspension. 



Bicerano et al. [30] examined the viscosity of suspensions of differ- 
ent hard bodies, and proposed a formula for the relative viscosity 
that provides a smooth transition between the dilute and the 
concentrated, Eq.([5|), regimes and is valid for low-shear 



V (0) = ^0 ( 1 
with 

and 



I in \ I in \ 

(6) 



1 + C^l - +C2 



Ci = M 4>c - 2, (7) 



C2 = kucfi - 2 M 0, + 1. (8) 

This theoretical formula accounts for both the low density exact 
results from hydrodynamics as well as the large volume fraction 
semi-empirical expansions. Pryamitsyn and Ganesan [13] pro- 
posed to extend the domain of validity of this expression to arbi- 
trary shear rates by considering the viscosity percolation thresh- 
old 0*(Pe) a function of the Peclet number. 

Recently, a model has been proposed that introduces in appro- 
priate form the correlations introduced by the excluded volume 
effects and gives an excellent quantitative description of the vis- 
cosity of solid and liquid suspensions of spherical particles at ar- 
bitrary filling fractions Ref. [5§M}] . This model incorporates an 
effective filling fraction 0e// that leads to an universal represen- 
tation of all experimental results in a master curve that suggest 
the role of 0e// as a scaling variable for the viscosity of these 
systems. 

Despite the importance of the theoretical and numerical analysis 
performed mostly in the case of spherical particles, in real systems 
the assumption of sphericity of the suspended particles is not al- 
ways satisfied since polydispersity and different particle shapes 
have to be taken into account [31]. Thus, the main objective of 



this article is to propose a continuum-medium description for the 



viscosity-concentration relation for different particle shapes as 
long as they are not too elongated. The description is derived us- 
ing DEMT techniques introducing correlations between particles 
through an effective volume fraction that incorporates excluded 
volume effects. These effects are responsible for the scaling prop- 
erties of the suspensions and show to be universal independently 
of the shear rate and the shape of the particles. 

The article is organized as follows, in Section II we propose a 
first correction to Eq.([T]) that takes into account excluded volume 
effects and use it as the starting point of a differential effective 
medium approach. In Section III we compare the predictions of 
our model with previous theories and with various experimental 
results for different particle shapes. Finally Section IV is devoted 
to conclusions. 



2 Correlations and DEMT approach 

The main difficulty when dealing with suspensions at large con- 
centrations is to take into account the correlations among par- 
ticles. Although many efforts have been devoted to microscopi- 
cally calculate the viscosity of a suspension of hard particles tak- 
ing into account the hydrodynamic interactions, the enormous 
mathematical complications associated to the many body prob- 
lem only permit to find corrections applicable to the low concen- 
tration regime. However, in first approximation, such correlations 
can be considered as follows: the contribution by the particles to 
the total stress tensor 11^ of the suspension is given in terms of 
an average over the volume V of the system in the form 

nr^^/n<'W, (9) 



P ~ y J — P 

where we represented the single particle contribution to the stress 
tensor by n^^-*. The upper V in Eq.(9) stands for the volume 
average and the factor A^ accounts for the contribution of the N 



independent particles. However, as defined in Eq.([9j), the average 
is strictly valid only for a system of point particles jSSl • 



Considering that each particle has a volume Vp, then the volume 
average must be performed over the free volume accessible to the 
particles defined by: Vjree = V — cNVp. Here c is a geometric fac- 
tor that takes into account the fact that the complete free volume 
cannot be filled with particles. Note that, for different shapes of 
the suspended particles, the value of the constant c will be dif- 
ferent. It also contains information about the maximum packing 
of particles the system may allocate and if the excluded volume 
effects are taken into account, the suspended phase contribution 
to the stress tensor is now given by [32], |iD] 






^v 



l-cc 



M0V^§ , (10) 



^v 



where = NVp/V, Vvq is the volume average of the traceless 
velocity gradient with vq the velocity field of the background fiuid 
f5D] . [55] . For finite-sized particles, this relation leads to the result 
that Einstein's expression scales with the excluded volume factor 
0/(1 — c0) instead of 0, and thus gives the following expression 
for the viscosity of a suspension 

r/(0) = r/o(l + W0e//), (11) 

where the effective filling fraction 0e// is defined by 



^eff 



1 — CC 



(12) 



The constant c depends on the filling fraction 0^ which is the 
critical concentration at which the suspension loses its fiuidity. 
This viscosity percolation threshold generally is greater than the 
purely geometrical percolation threshold of the particles, but less 
than or equal to the maximum packing fraction ^max [42] and is 



given by 

c=— - — . (13) 

0c 

The effective filling fraction (12) approaches the bare at low 



concentrations and becomes 1 at the divergence of the viscosity 
which occurs at 0c- The fact that the particles can not occupy 
all the volume of the sample due to geometrical restrictions is 
taken into account in the crowding factor c. For example, for a 
face centered cubic (FCC) arrangement of identical spheres, the 
maximum volume that the spheres may occupy is larger than for 
a random arrangement of spheres. 

The formula for the effective viscosity of the suspensions given 



by relation Eq. (|Tl|) can also be written as 

/ rh 

1 + M 



r] (0) = ^0 



,1 



(14) 



Although in this expression it is clear that 77 (0) incorporates the 
excluded volume corrections for the viscosity of a suspension, hy- 
drodynamic interactions are ignored and therefore one expects its 
validity been restricted to low concentrations. To improve it, fur- 
ther corrections must appear due to the interactions between par- 
ticles. In the system under consideration, these interactions are 
the hydrodynamic interactions which become increasingly impor- 
tant when increasing the filling fraction. The mentioned correla- 
tions can be accounted for by using DEMT techniques [12] . This 
theoretical method is based on a progressive addition of spheres 
to the sample in which the new particles interact in an effective 
way with those added in previous stages pi?] . 



Taking Eq. (11) as the starting point, suppose that we increase 
by (50e// the particle concentration in the suspension of viscosity 
rj {4>eff) by adding a small quantity ^4>eff of few new particles. 
If we treat the suspension into which we add these particles as 
a homogeneous effective medium of viscosity 77 (0e//), then the 



new viscosity can be written as 



T] {(j)eff + (50e//) = T] {(jy^f f) (1 + [j]] A0e//) . 



(15) 



Note that the increase in the effective particle concentration ^0e// 
is different from the effective concentration of new particles added 
at a given stage A0e//- This is due to the fact that one has to 
remove part of the effective medium, which already contains some 
particles, in order to allocate the new particles. From this, it 
follows that the fraction of particles of the new effective medium 
is given by (/)e// + ^0e// = 0e//(l - A0e//) + A(/)e//, from which 
we find 

^(Peff - ^ T • 

1 - y^eff 



(16) 



Substituting Eq.(16) in Eq.(15) and integrating we finally obtain 



r/(0) =r/o(l-0e//)" 



M 



(17) 



or, using the definition of 0e// 



1 (0) = ^0 



1 



,1 — C( 



-W 



(18) 



This relation for the effective viscosity of a suspension of rigid 
particles constitutes a powerful improved generalization of previ- 
ous theoretical results and empirical proposals, as we will show in 



the following section. In the limit of low concentrations Eq.(18) 



always reduces to Eq.([T) and the predicted Huggings coefficient, 
obtained by expanding Eq.(18) in a virial series, is given by 



/^ff = 2 '''1 ('''1 + 



(19) 



that has the form suggested by Douglas (see Eq.(17) of Ref . |12] ) . 



Although the procedure that leads to Eq. (18) is similar to other 



DEMT [31] -[37], [42], the introduction of excluded volume correla- 



tions through the effective filling fraction 0e// defined in Eq.(12) 
and then its role as integration variable in the differential proce- 
dure improves remarkably the agreement with experimental data 



when compared to other models, this is shown in Ref.jSS] for the 
special case of hard spheres. This agreement is due to the fact 
that the use of 0e// as integration variable implicitly considers 
correlations between spheres of the same recursive stage in con- 
trast to other models. 

It will be shown in the next sections that our model gives ex- 
cellent quantitative results at the whole concentration range for 
different particle shapes. However, at this point it is interesting to 



compare the functional form of our expression (18) with the one 



by Bicerano and coworkers, Eq.([6|). Their equation was designed 
explicitly to provide a smooth transition between the semidilute 



and concentrated regimes. Equation (18) can be written in a form 



similar to Bicerano's expression to obtain 



^0 



= 1 



1 + Ci 



+ C2 



+ C3 



+ 



(20) 

Here, 0c can be considered a function of the Peclet number and 
therefore is no restricted to the low-shear regime. The constants 
Ci and C2 are identical to the ones corresponding to the model 
of Bicerano and coworkers given by Eqs.([7]) and (|8]), and the 
coefficient of the cubic term is given by 



C3 = f W(M -3M + 2 



(21) 



It is not surprising that the coefficients of the linear and quadratic 
terms are the same in both models since by construction they re- 
duce to Eq.([I]) at low concentrations. Note however that in prac- 
tice, when using the model of Bicerano and coworkers, Eq.(|6]), 
once [77] is known for a given particle shape, 0c and kn are treated 
as two independent fitting parameters when comparing to exper- 
imental data 142]. In our model, on the other hand, 0c and kn 



are related by means of Eq.(19), therefore leading to one fitting 



parameter only. It is important to point out that although similar 

^c as 



in form, expression (20) has not the same divergence near 



6^ since the expression in square brackets in Eq.(20) diverges as 



.. However, if we make the expansion 



where 
and 



-(22) 
C[ = M ct>c - W\ , (23) 



C, = kHcfi - 2 M 0e (M - i) + M (M - 1) , (24) 

then, the expression in square brackets converges. This is an im- 
portant result meaning that the viscosity diverges near 0c with 
an scahng that depends on [rj] and therefore is not universal but 
depends on the particle shape and shear rate. 



3 Illustrative calculations 

3. 1 Spheres 

The first illustrative example of the procedure is a system of 
hard spheres which corresponds to take the value [rj] = 5/2 in 
the above expressions which leads to the well known Einstein's 
result, valid for very low concentrations 

V (0) = ^0 (l + 2^) . (25) 

Additionally to the models presented in the previous section, 
other phenomenological formulas have been proposed in order 
to fit experiments in the largest possible range of volume frac- 
tions. For example, Clercx and Schram obtained the following 
expression for the high-frequency effective viscosity [45J 



^oo (0) = ^0 



^^ 1^ + 1.4202 



1 - 1.420 
by considering two-particle hydrodynamic interactions only. 



(26) 
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Fig. 1. Left panels: (a) Relative viscosity i] (</>) /770 as predicted by our model 
[Eq.(18)] for spheres ([ry] = 5/2) at low- and high-shear rates as well as high- 
-frequency as a function of the volume fraction (/>. The lines correspond to the 
predictions of our model with RCP 4>c = 0.63 (upper line), close packing at FCC 
(\)c = 0.7404 (middle line), and fitting parameter (^^ = 0.8678 (lower line). The mea- 



sured data are from Refs. [l6] (circles) and [U] (triangles). Squares, Zhu et al. 
circles, van der Werff et al. [H]; triangles, Cichocki and Felderhof [20]. (b) Rep- 
resentation of various viscosity data as suggested by Bedeaux (Ref. [23]). Squares, 
SJ18 low-shear limit (Refs. [1] and [IS]); circles, SJ18 high-shear limit (Refs. [26] 
and |46j ) : triangles, high-frequency limit of the real part of the complex shear vis- 
cosity (Ref.[38j). The lines are the results of our model with (pc = 0.63 (upper line), 
(j)c = 0.7404 (medium line), and (pc = 0.8678 (lower line). Right panels: The same 
as the left panels but the lines correspond to the model by Krieger and Dougherty 
[Eq.li)]. 

In order to carry out comparisons of our model with experiments 
and other models it is important to notice that the value of 0c in 



Eq.(13) is a free parameter of the theory to be chosen in order 



to best fit the experimental results. Nonetheless, this parameter 
can be chosen beforehand based on physical arguments, and then 
used to compare with specific experimental situations. 



This is done in the left panel of Fig. 1, where the behavior of the 
relative viscosity 77 (0) /t]q (Fig. la, left) is compared with exper- 
imental results of de Kruif et al. [IB] and of Krieger [ITj for low 
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and high-shear rates, and with the results of van der Werff et at. 
p[5] and by Zhu et at. [IS] at high frequencies as a function of the 
volume fraction 0. We also plot the values obtained by Cichocki 
and Felderhof |[2D] for the high-frequency case. The upper curve 



represents the prediction of Eq.(18) with 0^ = 0.637, which cor- 



responds to the random close packing (RCP) of identical spheres. 
The comparison with the experimental results for low-shear rates 



is excellent. The middle curve is the prediction of Eq.(18) with 
(j)c = 0.7404, which corresponds to FCC close packing. This gives 
again an excellent agreement with the experimental results for 
the case of high-shear rates. These results are consistent with 
the known fact that for low-shear rates the spheres remain dis- 
ordered while at high-shear rates the equilibrium microstructure 
of the dispersion is completely destroyed and the spheres adopt 
an ordered FCC configuration. The lower curve is the prediction 



of Eq.(18) with 0^ = 0.8678, which gives again an excellent fit 
with the infinite-frequency data and with the values obtained by 
Cichocki and Felderhof 



In order to account for thermodynamic interactions between the 
spheres, Bedeaux proposed the following expression for the vis- 
cosity [23]-[21] 

mhi_l = Hi + sm, (27) 

^ (0) /^O + 2 

where *5'(0) is an unknown function of the volume fraction giv- 
ing the modification of the moment of the friction forces on the 
surface of a single sphere due to the ensemble-averaged hydrody- 
namic interactions with the other spheres [ [24] . 

As discussed by Bedeaux [23], S (0) is a more sensitive represen- 
tation of the relative viscosity data since its expansion in powers 
of (j) converges much better. For this reason, the differences be- 
tween the predictions of the different models are more noticeable 
in this representation. In Fig. lb, left panel, we plot the function 
*S'(0) for various experimental results obtained in a variety of ex- 
perimental situations and compare with the values obtained with 
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Fig. 2. Left panels: The same as in Fig.l but the hnes correspond to the model by 
Quemada [Eq.([5])]. Right panels: The same as in Fig. 1 but the lines correspond to 
the model by Bicerano and collaborators [Eq.Q]. We also show the high-frequency 



result by Clercx [Eq.(26)] 
(dash-dotted line). 



our model. We take the same values of 0c used previously for the 
low-shear, high-shear, and high-frequency cases. The agreement 
with the experimental data is very good especially at large values 
of 0. Note that for small values of the experimental accuracy 
of 5* (0) is unsatisfactory which explains the large scatter in the 
experimental points 



As a comparative, in the right panel of Fig. 1 and in Fig. 2 we 
show the predictions obtained with the models of Krieger and 
Dougherty [Eq.([3)], Quemada [Eq.([5])], Bicerano [Eq.(|6])], and 
Clercx [Eq.(26])]. In each case, the symbols (a) and (b) state for 
the direct comparison with experiments (a) and after using Be- 
deaux's expression (b). It is clear that our model matches the 
data at all volume fractions better than the other models tested. 

Let us stress here that the improvement of our model given 
by Eq.(18) to the models given by Eqs.([3]) and ^ is due to 
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the use of 0e// given by Eq.(12). Indeed, the model of Krieger 



and Dougherty can be obtained from the dilute-Umit expression 
ij (0) = r/o (1 + [r/] (f)) by introducing the effective filhng fraction 
0^y^ = 0/0max, which is larger than 0. This definition of 0^^^ un- 
derestimates the available volume for the particles at low (and 
therefore, overestimates 0e//) while tends to the correct limit at 
high (J). In order to obtain the correct dilute limit, the overesti- 
mation of 0^/, has to be compensated by decreasing the hydro- 
dynamic drag factor by the same constant factor 0niax, that is, 
[77] = [77] 0max- Then, the dilute-limit expression can be rewrit- 
ten as 

Vr{<l>) = (l + lvf"'Pff'j)- (28) 



Now, Eq. (|3]) can be derived from (28) by following the differen- 
tial method used in the previous section with the effective filling 
fraction 0^j instead of 0e// and the corresponding [77] instead 
of [77]. The relative difference {4>^/f — 4'eff)/4'eff is a decreasing 
function of that vanishes at 0max (we assumed that 0c = 0max 
in this analysis. Therefore, the overestimation of the filling frac- 
tion in 4>fff, is progressively less important with increasing and 
the constant underestimation of the hydrodynamic drag term in 
Krieger and Dougherty's model cannot be compensated by (j)f/f- 
Thus, Krieger and Dougherty's model underestimates the viscos- 
ity of the suspension at large volume fractions, as confirmed in 
Fig. lb. The same reasoning can be applied to Bullard's model 
since 1/K = (jy^ plays the role of a critical concentration and there 
is no mathematical distinction between Eq.(|4]) and the Krieger 
and Dougherty equation. 

In what follows, we show a direct comparison of the functional 
forms of our proposal and the other models as well as a virial 
expansion of them. For example, at high-frequencies, a virial ex- 



pansion of Clercx and Schram's expression, Eq.(|26[) gives 

'- ^ = l + -( 

1 - 1.420 2 



r/00 (0) /m = l+ 'r .\Z = 1 + ^,0 + 4.9702 + 7.0603 + 



(29) 
while our expression, Eq.([l8|) with [7/] = 5/2 and 0c = 0.8678, 
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can be expressed using a Pade approximation as 

-700 {^) /% = i+ i-^^ + l-y = 1 + ^0+4.76^^+7.950^+0 (0*) , 

(30) 
showing that the second and third virial coefficients are very 
close. Similarly, the expression by Bicerano and coworkers, Eq.([6]), 
can be written in the case or hard-spheres with 0c = 0.637 as 

1 (0) /m = (l - ^) [1 - 0.6280 + 0.840^; 

= 1 + ^0 + 6.2602 + 13.4703 + O (04) , (31) 
while our expression gives for the same value of 0c 



-2 



1 (0) /m =[l- ^ j [1 - 0.640 + 0.4150 

= 1 + ^0 + 5.802 ^ 12.360^ + O (04) . (32) 

Both expressions are very similar up to third order in 0. However, 
as explained before, near 0c both models predict different asymp- 
totic relations. While Bicerano's model predicts a divergence of 
the viscosity as (1 — 0/0c) , our model predicts that the viscos- 



ity diverges as (1 — 0/0c)^ ' , Eq.(22), for the case of spherical 
particles. Our prediction is in close agreement with that of Cheng 
et. al.jlil, that adapted mode coupling theories developed for su- 



percooled liquids and the liquid-glass transition for molecular sys- 
tems to calculate the low-shear viscosity of colloidal dispersions. 
Taking into account hydrodynamic interactions, they found that 
near the glass transition 0^ = 0.62, the viscosity diverges approx- 
imately as (1 — (j)/(j)g)~ ' which is in close agreement with our 
model if 0c — (j)g. 

The proposal by Pryamitsyn and Ganesan[j43j to extend Bicer- 
ano's expression to arbitrary shear rates by considering the pack- 
ing fraction 0c(Pe) a function of the Peclet number can also be 
tested by expanding in series of for a large shear rate, which 
means to take 0c = 0.7404. In this case Pryamitsyn and Ganesan 
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Fig. 3. Left panel: Intrinsic viscosity for ellipsoids as a function of the aspect ratio 
p = b/a. Right panel: Relative viscosity rj (cj)) /r]o as predicted by our model for 
ellipsoids with different values of p. The critical packings were assumed to coincide 
with the MRJ states as obtained in Refs. ([52j) and (L53j): for p = 1/4, (j)c ^ 0.64; 
for p = 1/2, (pc — 0.7; for p = 1, 4 
(hr ~ 0.632. 



0.637; for p = 2, ^c - 0.704; and for p = 4, 



obtain 



V (0) /m 



1 



-2 



1-0.540 + 0.6220^ 



0.7404; 
= 1 + 1.60 + 3.102 ^ 5_4g^3 ^ Q f^4 



(33) 



while our model gives 



1 (0) /^o = 1 



0.7404. 



1- 0.2010 + 0.3220^ 



= 1 + -0 + 5.250^ + 9.940^ + O 



(34) 



Notice that the agreement is not very good even at first order 
in where the extended Bicerano's expression do not agree with 
Einstein's expression for low 0. 
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L/a=1,b/a=1 
L/a=0.74, b/a=1 .26 
Lya=0.5, b/a=0.5 
Lya=0.2, b/a=0.£ 




Fig. 4. Left panel: Relative viscosity ij {(p) /rjo as predicted by our model for 
flat-ended cylinders with different aspect ratios. Right panel: The same as left panel 
but for hemiellipsoidal-ended cylinders with different aspect ratios. In both cases 
we assumed that the critical packing coincides with the random packing of sphe- 
rocylinders with the same aspect ratio [ 



3.2 Ellipsoids 



Another case in which the behavior of highly concentrated sus- 
pensions is of importance corresponds to suspensions of elhp- 
soidal particles, since they have been widely used as models of 
globular proteins and also because the hydrodynamic properties 
of ellipsoids are analytically known. For example, the intrinsic 
viscosity of an ellipsoid is given by [50] 



[v] = Up'-i 



(V _ i) 5 + 2j92 + 1 



+ 



3p^ {3B + 2p^ - 5) [(2j92 + 1)5-3] 

14 2 

+ 



3p2 (35 + 2p^ - 5) {p^ + 1) {-3p^B + j?2 + 2) 
p'^-l 1 

^ p^ {p^ + 1) [(2j92 _ 1) 5 _ 1] I ' (^^) 
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L/a 


1/4 


1/2 


1 


2 


4 


[v] 


3.92 


3.05 


2.87 


3.40 


5.27 


M^) 


0.67 


0.69 


0.68 


0.62 


0.53 



Table 1 

Intrinsic viscosity [rj] [taken from Ref. ([8J)] and critical packing cpc [taken from Ref. 
([5^)] of flat-ended circular cylinders. 

where 



,-1 ^^2 



B = p Wj9^ — l) cosh j9, when 



= 1, 

— 1 /i 2\~l/2 —1 

= P [i- — P ) COS p, 



P>1, 

when J3 = 1, 

when p < 1. (36) 



Here, the axis ratio is defined by p = b/a with b the polar ra- 
dius and a is the equatorial radius. In the left panel of Fig. 3 
we show the behavior of the intrinsic viscosity as a function of 
the axis ratio p. As can be seen, the lower value of the intrin- 
sic viscosity occurs for spherical particles and increases slightly 
more sharply for prolate ellipsoids than for oblate ones. In the 
right panel of Fig. 3 we represent the viscosity as a function of 
the concentration for ellipsoids with different values of p. In all 
cases we have assumed that the critical packing coincides with 
the so called maximally random jammed (MR J) state which cor- 
responds to the least ordered among all jammed packings ^51j . 
The critical densities of simulated packings of elfipsoids are cal- 
culated in Refs. ([52]) and ([l53j) for different aspect ratios and 
we used them in the results of Fig. 3. As can be seen, the non- 
monotonic behavior of the intrinsic viscosity shown in the left 
panel of Fig. 3 is reflected in the right panel of this flgure. 

3. 3 Cylinders 



Cylinders may be used to represent DNA molecules in certain 
physical conditions and also as models for rod-shaped viruses. 
Flat cylinders are encountered most often. Except in certain limit, 
analytical results are not available for hydrodynamic properties 



L/a 


1 


0.74 


0.5 


0.2 


b/a 


1 


1.26 


0.5 


0.8 


H 


2.944 


2.910 


2.571 


2.512 


</.,(«) 


0.62 


0.62 


0.68 


0.68 



Table 2 

Intrinsic viscosity [r/] [taken from Ref. ([8])] and critical packing 

([53])] of hemiellipsoidal-ended circular cylinders. 



[taken from Ref. 



of cylinders. Thus, numerical calculations have to be used. In the 
limit that the total length (2L) goes to zero, a flat-ended cyhnder 
becomes a disk. The intrinsic viscosity of a disk with radius a can 
be found from Eq.(35) by taking the j9 — )► limit. The result is 
[j]] = 128a^/45. In the case of globular cylinders for which the as- 
pect ratio L/a ranges from 1/4 to 4 the hydrodynamic properties 
obtained numerically are listed in Table 1, [jSj. Together with the 
hydrodynamic properties, in the last column of Table 1, we show 
the critical packing which we assume to coincide with the ran- 
dom packing of the cylinders. Actually, we are not aware of data 
for the close packing of flat-ended cylinders so that we are using 
the results for spherocylinders obtained in Ref. (fS^) with the 
same aspect ratio. The viscosity as a function of concentration is 
shown in Fig. 4. 

In the case of circular cylinders with hemiellipsoidal ends, the 
results for the intrinsic viscosity jH] and the critical packing are 
listed in Table 2. As in the previous case, the value of the random 
packing listed corresponds to spherocylinders with the same as- 
pect ratio defined as (6 + L) /a (see inset of Fig. 5). The viscosity 
in this case is shown in Fig. 5 where one can see that it is larger 
for the cylinders with the largest aspect ratio. 

Experimental data for the concentration and particle size-dependence 
of the low-shear viscosity of isotropic rod-dispersions are dis- 
cussed in Ref. [9]. Four systems of stiff rods are considered: col- 
loidal silicaboehmite, xanthan (A = 0.5), schizophyllan (A = 0.3) 
and non-Brownian PMMA-fibre. Here, A is the ratio of the poly- 
mer length and the persistence length. The measured intrinsic vis- 
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Silica rods 


Schizophylian 


Xanthan 


PMMA-fibre 


iv] 


50.2 


42.7 


58.1 


27.6 


L/a 


44.0 


46.0 


56.0 


39.8 


M^) 


0.245 


0.235 


0.193 


0.271 



Table 3 

Table 3. Intrinsic viscosity [ij] [taken from Ref. (|9])], aspect ratio L/a, and critical 

packing 



obtained from Eq.(37)] for isotropic dispersions of rods. 



cosities [77] of these systems are tabulated in Table 3. The packing 
density is determined by their aspect ratio. For a random packing 
of thin hard rods, it was shown[^,|[5S] that 



^max- ='^ /or - >> 1, 



a 



a 



(37) 



where ^ is the average number of contacts experienced by a rod. 
Experiments |E5], [STj on random rod packing yield ^/2 = 5.4. 
This isotropic maximum packing fraction is a metastable glass 
with respect to the thermodynamically more favorable nematic 
phase as predicted by Onsager. The aspect ratio and the cor- 



responding maximum packing calculated using Eq.(37) for the 



systems considered in Ref. P, are also tabulated in Table 3. In 
Fig. 5a we compare these data to our model using as parame- 
ters the values given in Table 3. As can be seen, poor agreement 
is found between the model and the experiment. A number of 
reasons can be argued to explain this result. For example, the 
silica rods are weakly attractive which implies a steeper viscosity- 
concentration curve as compared to the rigid macromolecules 
xanthan and schizophylian. Of these two rod-like macromolecules 
the xanthan chain possesses higher flexibility which may also in- 
fluence the viscosity data. Additionally, since long rod-shaped 
molecules "entangle" in the semidilute regime, it is expected that 
some orientational correlation arises at large concentrations and 
our model do not consider these effects. 

To explore the possibility that the rods present some amount of 
orientational correlations, we have compared our model to the ex- 
perimental data but considering that some alignment due to the 
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Fig. 5. Left panel: Relative viscosity r\ {(p) /r]o as predicted by our model for isotropic 
rods with different aspect ratios and intrinsic viscosities as given by Table 3. Silica 
rods (solid line), schizophylian (dashed line), xanthan (dotted line), and PMMA 
fibre (dash-dotted line). Right panel: The same as left panel but assuming orienta- 
tional order that leads to 0max — 0.9069. The fitted values of [q] were 90 for silica 
rods (solid line), 50 for schizophylian (dashed line), 40 for xanthan (dotted line), 
and 13 for PMMA fibre (dash-dotted line). The data sets were taken from Ref.^9j. 



flow and the correlations is possible. In this case, the maximum 



packing is no longer given by Eq.(37) and we assume the oppo- 



site limit, that is, that the rods are completely ahgned. Then the 
maximum packing takes the value corresponding to an hexagonal 
arrangement of parallel rods 0max — 0.9069. Then, the intrinsic 
viscosities for random rods as tabulated in Table 3 are no longer 
useful and we take [77] as a fitting parameter. This comparison is 
shown in in Fig. 5b. As can be seen, a much better agreement 
is obtained using these assumptions. Thus, some degree of align- 
ment is suggested by the model. Nonetheless, for a concentration 
dependent orientational correlation, the intrinsic viscosity [77] is 
rather a function of and therefore, the values predicted by the 
present model may not be correct. 
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3-4 Dumbbells 

A dumbbell that consists of two spheres can be used as a model 
for a protein dimer or a protein consisting of two separate do- 
mains. Wakiya [7] and Brenner [B] calculated the intrinsic viscos- 
ity for dumbbells consisting of two equal-radius spheres at various 
separations. Here, we will only discuss the case when the ratio 
L/a = 1, where a is the radius of the spheres and L is half of 
their center-to-center separation. In this case, [r/] = 3.4496. 

The viscosity of dumbbells made of two fused spheres is shown in 
Fig. 6 and compared with the numerical results reported by in't 
Veld and coworkers for nanodimers [i58j . The closed circles repre- 
sent points assuming an effective volume fraction for particle radii 
adjusted to the peak onset in the nanoparticle pair distribution 
function. This adjustment corrects for the solvation shell around 
the nanoparticles when the solvent is treated explicitly [5B] . The 
fitting parameter in this case was 0c — 0.58, which is close to the 
critical packing value for the glass transition predicted by mode 
coupling theory 0c — 0.56 (see Ref. [52]). 



3.5 Other Shapes 

Notice that the intrinsic viscosity [r/] for long rod-shaped particles 
takes values much larger as compared to the corresponding ones 
for spheres. However, it is possible to have very large values of 
the intrinsic viscosity without making a very extended or flat 
object [32]. An strategy to increase [77] is to consider irregularly 
shaped particles like sponges or jack-like objects. The intrinsic 
viscosities for these shapes have been calculated numerically by 
finite element computations in Ref. [32]. We use these values in 
Fig. 7 to compare the viscosity-concentration curves for three 
representative cases: a sponge, a wire frame, a square ring, and 
a jack-like object. 
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0.0 0.1 0.2 0.3 



^ 



0.4 0.5 



Fig. 6. Relative viscosity rj {(p) /r]Q as predicted by our model for a system of dumb- 
bells. The closed squares represent points assuming radii of the spheres a = 5a, 
where a is the size of a Lennard-Jones solvent atom in determining (j), whereas 
the circles denote a volume fraction for radii adjusted to the peak onset in the 
nanoparticle pair distribution function. 



The sponge is constructed starting with a cube in which a square 
channel is cut through the center of each face, which passes com- 
pletely through the cube, as seen in Fig. 7. The parameter m 
is taken to be the edge length of the cutout face in units of 
the cube edge length. When m approaches 1 a rigid cubic wire 
frame is obtained. A similar procedure is employed to construct 
the flat square. The jack is constructed by poking three rectan- 
gular parallelepipeds orthogonally through a sphere. In all the 
curves of Fig. 7 we assumed 0c — 0.637, the random close pack- 
ing for spheres and the values of the intrinsic viscosity given in 
Ref.pS]. As expected, the objects with larger intrinsic viscosities 
have steeper curves. These examples are intended for illustrative 
purposes only, since we are unaware of experimental results to 
compare with. 
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4 Conclusions 



We have presented a simple model based on an effective- medium 
theory for the calculation of the viscosity of suspensions of arbitrarily- 
shaped particles as a function of particle concentration. The model 
considers excluded volume interactions between the particles through 
an effective filling fraction (peff- This quantity introduces a uni- 
versal scaling that may be used to reduce both experimental and 
theoretical results to a master curve |22],|1D] which is indepen- 
dent of the experimental details or the shape of the particles. 

Starting from known values and formulas for the intrinsic viscos- 
ity of the particles, the procedure yields to analytical expressions 
that predict the viscosity of the system for the whole range of con- 
centrations. At low filling fractions it reduces to the correct limit 
while at high concentrations it diverges in a way similar to that 
predicted by mode coupling theories. In contrast to other mod- 
els jBO], our proposal contains only one fitting parameter which 
corresponds to the critical packing where the suspension loses its 
fluidity. 

When applied to a suspension of spherical particles, our model 
improves considerably the predictions obtained using the well 
known Krieger and Dougherty model and any other model tested 
in the whole concentration range. We have employed our model to 
predict the viscosity of elhptical, and cyhndrical particles, as well 
as dumbbells made of fused spheres and other complex shapes. 
In all cases where numerical or experimental data are available, 
the agreement with the proposed model is very good. It is con- 
venient to emphasize that our model is not intended to describe 
correctly suspensions of large fibres, since in this case orienta- 
tional correlations may exist. These correlations could introduce 
dependences of the intrinsic viscosities on the concentration. In 
a previous work [ID] we have applied the procedure to emulsions 
of spherical droplets with equally good results. 
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Fig. 7. Relative viscosity r/ {(p) /rjQ as predicted by our model for a sponge 
(m = 23/27, [r]] = 44.7), a cubic wire frame (m = 33/35, [rj] = 255) a square ring 
(m = 23/25, [rj] = 98.7), and a jack-like object ([r/] = 3.68). The specific propor- 
tions of the jack are as follows, if the width of the parallelipipeds have unit length, 
the length has 15 units and the sphere diameter has 9. In all cases we assumed 
(t)c ^ 0.637. 



Due to the importance of shape effects on the rheological behav- 
ior of colloidal dispersions and despite that there are numerous 
reports for industrial systems but fewer data for suspensions with 
controlled geometry, we consider that the results presented in this 
article can help to provide a valuable characterization of these 
systems with very promising practical applications. 
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